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ABSTRACT 

We report the most rapid rate of period change measured to date for a pulsating DA (hydrogen 
atmosphere) white dwarf (WD), observed in the 292.9 s mode of WD 0111+0018. The observed 
period change, faster than 10~^^ s s~^, exceeds by more than two orders of magnitude the expected 
rate from coohng alone for this class of slow and simply evolving pulsating WDs. This result indicates 
the presence of an additional timescale for period evolution in these pulsating objects. We also 
measure the rates of period change of nonlinear combination frequencies and show that they share 
the evolutionary characteristics of their parent modes, confirming that these combination frequencies 
are not independent modes but rather artifacts of some nonlinear distortion in the outer layers of the 
star. 

Subject headings: stars: white dwarfs-stars: individual (WD 0111+0018)-stars: oscillations (including 
pulsations)-stars: variables: general 
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1. INTRODUCTION 

White dwarf (WD) stars represent the final evolution- 
ary stage of all single low-mass stars, and are thus rep- 
resentative of the future of the majority of stars in our 
Galaxy. They are primarily composed of the inert by- 
products of hydrogen and helium fusion, and their evo- 
lution is dictated by the rate at which these carb on and 
oxygen ions lose their residual thermal energy (|Meste]| 
[l952). This cooling takes billions of years. 

However, we can witness this cooling on human 
timescales by watching the WDs that pulsate. These 
variable stars come in several classes, based on their 
dominant atmospheric composition (their strong surface 
gravity and rapid gravitational settling leads to chemi- 
cally homogenous photospheres). Here we will restrict 
our discussion to the coolest class of such pulsators: the 
hydrogen-atmosphere DAV (or ZZ Ceti) stars, found in 
an instability strip betwe en 10800 — 12300 K for standard 
log q = 8.0 (IKoester Holbcrg 2001; Bergero n et al.l 



2664'; 'Mukadam et al."20M; iGianninas et al.ii2005n200e 
Mukadam et al. 2006). 

Aside from their variability, which is brought on by 
a hydrogen partial-ionization zone in the star's non- 
degenerate atmosphere, DAVs appear to be normal WDs. 
They are therefore belie ved to be r epresenta tive of evo- 
lution for all DA W Ds (iRobinsonI [l979: Fon taine et"all 
[T985H Bct geron et al.ll2004[ ) . Cooling from neutrino emis- 
sion is expected to be ne gligible for DA WDs within the 
ZZ Ceti instability strip (jWinget et al.ll20M ). and their 
evolution should be dominated by radiative surface emis- 
sion. 

The C/O-core DAVs undergo multi-periodic, non- 
radial g-mode pulsations with periods between roughly 
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100 — 1400 s. Seismology using these g-mode pul- 
sations has enabled us to constrain the mass, core 
and envelope composition, rotation rate, and the be- 
havior of convection in these objects (see review s by 
Winget & KeDlei< ^2001 iFontaine fc BrassardI [20081 and 
Althaus et al. 2010!). These modes are often stable 
in period and amplitude, especially for the hot DAVs 
(hDAVs) near the blue-edge of the ZZ Ceti instability 
strip. hDAVs are expected to show an extremely slow 
period drift caused by the gradual cooling of the star. 

Two of the longest-studied hDAVs are G117-B15A and 
ZZ Ceti itself (also referred to as R548), which have been 
observed since the early 1970s. Both show considerable 
stability in their largest-amplitude modes. The rate of 
change of period with time, dP/dt, of the dominant 215.2 
s mode in G117-B15A has been measured to be (4.19 ± 
0.73) X 10"i5 s s-i (jKepler et al.|[20Tll) . R548 has been 
observed since 1970, but more sparingly: Fo r the 213.1 s 
periodicity in R548. [Mukadam et al.l (|2012l ) determine a 
dP/rft = (3.3+1.1) X 10-15 ss-^ These values have taken 
decades of observations, and are in line with theoretical 
predictions of rates of period change fro m cooling alone 
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in DAVs of dP/dt < IQ-^-^ s s'^ (Bradlev et al 
IBischoff-Kim et al.l[2008l : ICorsico et al.ll2012D . 

After more than nine years of monitoring, we have 
found evidence for an hDAV that has a rate of pe- 
riod change inconsistent with cooling alone. The period 
change in the highest-amplitude periodicity for that star, 
SDSS JOlllOO.63+001807.2 (hereafter WD 0111+0018; 
3 = 18.7 mag), exceeds 10"^^ s s-^. This observed rate 
of change of period with time is several orders of mag- 
nitude faster than predicted from cooling alone for this 
11810 K WD. It therefore signifies some physical process 
operating on a considerably shorter timescale. 

In this paper we present our observations of period evo- 
lution in the four highest-amplitude periodicities present 
in WD 0111+0018. In Sections 2 and 3 we outline our 
observations and analysis, and we reserve Sections 4 and 
5 for a discussion of possible explanations for this rapid 
period change and our conclusions. 
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TABLE 1 
Journal of observations. 



Subgroup 


UT Date 


Length 


Seeing 


Exp. 






(hr) 


(") 


(s) 


2003a 


n Jan zOUo 


1.5 


2.3 


15 




Uz reb ZUUvS 


0.9 


1.6 


20 


2004a 


iNov zOUo 


2.9 


2.3 


15 




oU iNOV zUUo 


3.0 


2.4 


15 


2005a 


iZ Dec ZUU4 


5.0 


2.5 


15 




14 Dec zUU4 


4.2 


2.7 


15 




Id JJec ZUU4 


3.4 


1.4 


15 


200Da 


Uo Dec zUUo 


3.5 


1.8 


15 




Ud Uec zUOo 


5.4 


1.5 


15 


2006b 


oU Uec zOUo 


3.4 


2.0 


15 




ol Uec zUOo 


2.8 


2.0 


15 




Uo Jan zUUo 


2.6 


1.6 


15 


2007a 


ly bep zUUd 


3.1 


1.5 


15 




Z2 bep ZUUo 


2.7 


1.7 


15 


2007b 


ZO iNOV ZUUo 


6.1 


1.3 


10 


2010a 


16 Aug zuuy 


2.2 


1.6 


15 




z4 Aug zuuy 


4.2 


1.4 


15 


2010b 


15 Dcp zUU9 


5.5 


1.5 


15 


2010c 


11 iNov zuuy 


3.5 


1.6 


15 




i / iN ov zuuy 


o. / 


2.7 


30 




o T M-->, T onnn 

zi iNOV zuuy 


7.4 


1.6 


30 


2010d 


1 o T 1 -. oni n 
iz Jan ZUiU 


3.3 


2.1 


15 




1 Q T „ om n 
io Jan ZUiU 


3.6 


1.8 


15 


2011a 


1 o o„», om n 

\1 bep zUlU 


3.6 


1.1 


15 




T *? o om n 
Id bep ZUlU 


6.0 


1.6 


15 


2011b 


nv r^^i-*- om n 
U; Uct ZUiU 


4.5 


1.3 


15 




no f^^i- om n 
Uo Uct ZUiU 


4.3 


1.6 


15 




09 Oct 2010 


5.3 


1.3 


15 


2011c 


10 Nov 2010 


1.7 


1.8 


15 




11 Nov 2010 


3.1 


1.3 


15 


2011d 


02 Jan 2011 


2.8 


3.0 


20 




04 Jan 2011 


4.1 


1.6 


15 




05 Jan 2011 


1.7 


1.7 


20 


2012a 


02 Oct 2011 


3.8 


1.2 


15 




03 Oct 2011 


2.6 


1.2 


15 


2012b 


23 Oct 2011 


4.5 


1.5 


15 


2012c 


01 Feb 2012 


1.2 


1.4 


15 




02 Feb 2012 


1.2 


2.5 


15 




03 Feb 2012 


1.8 


2.1 


15 



2. TIME-SERIES PHOTOMETRIC OBSERVATIONS 

iMukadam et al.l (|2004[ ) discovered pulsations in 
WD 0111+0018, resulting from a search for variable 
DAs. The WD was initially targeted as a result of 
model atmosphere fits to spectra from the Sloan Digital 
Sky Survey: with T^s = 11510 ± 110 K and log^ = 
8.26 ± 0.06, the object was predicted and confirmed to 
be a DAY within the ZZ Ceti instability strip. 

Using new treatm ents of line broadening theory, 
iTremblav et al.l (|2011l ) refined the temperature and sur- 
face gravity of WD Olll-f-0018 to T^s = 11810 ± 190 K 
and log^ = 8.17 ± 0.07, which corresponds to a mass of 
MwD = 0.71 ± 0.04 Mq and a cooling age of 510 Myr. 

Including the initial discovery observations in early 
2003, the White Dwarf Group at the University of Texas 
at Austin has logged a total of more than 136 hr of obser- 
vations of WD 0111-1-0018 on 39 separate nights. As the 
journal of observations in Table[T]indicates, the data span 
more than nine years. One motivation for these observa- 
tions was to use the stable pulsation modes of hDAVs as 
precise clocks, searching for evidence of planetary com- 
panions by lo oking for periodic variations in the puls e 
arrival times (jWinget et al.l 120031: iMullallv et alll2008[ ). 
The project has also effectively opened a window on new 
evolutionary timescales for more than a dozen hDAVs, 
including WD 0111-^0018. 
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Fig. 1. — The pulsation spectrum of WD 01114-0018 for a typieal 
run, this on 15 Sep 2009. A Fourier transform of this 5.5-hr run 
illustrates that the four periodicities of interest, seen marked, can 
be resolved in a single night. 

Ah of our data on WD 0111+0018, more than 33,100 
images, have been taken using t he same instrument 
(Argo s, a frame-transfer CCD; see INather fc MukadamI 
120041 ) with the same filter (a 2 mm BG40, to reduce sky 
noise) on the same telescope (the 2.1m Otto Struve at 
McDonald Observatory), allowing for an especially co- 
herent data set. Depending on the conditions, we have 
used 10 — 30 s exposures, as indicated in Table [1] 

The raw science frames are calibrated by dark sub- 
traction and flat-fielding. We perform weighted, cir- 
cular, aperture photometry on the calibrated frames 
using the external IRAF package ccdJisp written by 
Antonio Kanaan (t he re duction method is o utlined in 
iKanaan et all 12003 and IMullallv et all 12001 . We di- 
vided the sky-subtracted light curves using two bright 
comparison stars in the field, SDSS J011055.57+001850.6 
{g = 15.0 mag) and SDSS J011056. 17+001955.8 [g = 
16.2 mag), to correct for t ransparency variations. Usin g 
the WQED software suite (jThompson fc Mullallvir2009D . 
we fit and subtract out a low-order polynomial (at the 
timescale of several hours) to each run to remove any 
long-term trend caused by atmospheric extinction, and 
apply a timing correction to each observation to account 
for the motio n of the Earth around the barycenter of the 
solar system (jStumpfllll980l) . 

3. LIGHT CURVE ANALYSIS 

3.1. The Pulsation Spectrum of WD 0111+0018 

The light curve of WD 0111+0018 features modulation 
at four distinct periodicities between 100 — 300 s, which 
can be well-resolved in a 3-hr run in good conditions (see 
Figured]). In order to most accurately identify these four 
periods, we begin by computing a Fourier transform of 
our entire data set. Using these results as an initial guess, 
we then perform a nonlinear least-squares fit for the fre- 
quency, amplitude, and phase of the four periodicities 
present. 

This allows us to isolate an initial estimate of the pe- 
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riods of the four highest-amphtude pulsations present in 
WD 0111+0018, which we use to construct an initial 
(O — C) diagram. The amplitudes for these four peri- 
odicities never fall below 0.5%; the amplitudes listed in 
Table [2] are weighted-mean amplitudes using each sub- 
group. 

We perform this analysis independently for all 18 sub- 
groups of our data. The subgroups are created by com- 
bining data from as many contiguous nights as possi- 
ble; no combined subgroup spans more than two weeks. 
These subgroups have been identified in Table [TJ In each 
case, the periods we determine for each subgroup match 
the periods we have found for the entire data set, within 
the uncertainties. 

We fix these periods and compute a simultaneous linear 
least-squares fit for each subgroup to create an (O — C) 
diagram. Thus, the uncertainties in our (O—C) diagrams 
represent formal uncertainties on the least-squares fit for 
the phase on at least one — and as many as three — nights 
of data, separated by as many as 10 days. 

We will refer to the two highest-amplitude periodici- 
ties as the parent modes /i and /2, and the two smaller- 
amplitude periodicities 2/i and /i + /2 as the combi- 
nation frequencies, such that the frequency for 2/i is 
(within the uncertainty) exactly twice the frequency of 
/i and the frequency for /i + /2 is the sum of the fre- 
quencies of /i and /2. We will discuss in more detail 
these combination frequencies and how they arise in Sec- 
tion H+l 

Additionally, we note the presence of an additional 
combination frequency 2/2 at 127.8 s, visible but not 
labeled in Figure [TJ In fact, we detect at least five 
additional nonlinear combination frequencies. Using 
all of our light curves from 2009 to early 2012 to de- 
tect as many combination frequencies as possible, we 
find: /i = 2.739(29)%, /2 = 1.526(31)%, 2/i = 
0.472(27)%, /i + /2 = 0.732(28)%, 2/2 = 0.193(25)%, 
3/1 = 0.143(44)%, 2/1 + /2 = 0.295(98)%, /i + 2/2 = 
0.147(44)%, 3/2 = 0.068(28)%. However, the amplitudes 
of most of the combination frequencies are too low to pro- 
duce an (O — C) diagram for all subgroups, so we do not 
include them in the (O — C) analysis. 

3.2. Constructing an (O — C) Diagram 

We demonstrate a secular change in the pulsation peri- 
ods of WD 0111+0018 by constructing (O — C) diagrams, 
where we compare the observed time of maximum for a 
pulsation (O) to when we expect such a maximum as- 
suming that the pulsation obey s a constant period (C). 

Following iKepler et al.l ()1991[ ). if the pulsation period 
is changing slowly with time, we can expand the observed 
time of maximum of the E*'' pulse, i^, in a Taylor series 
around Eq: 

tmax \e = tmax \eo H 77^ l-Eo ^ ^o) + 

ah 

1%!^ {E ~ Eo)' + ... (1) 

where the epoch is E = t/P and the change in arrival 
time with epoch, dt/dE, is the period, P. If we drop 
all terms higher than second order (assuming that P is 



negligible), we arrive at the classic (O — C) equation 

O - C = Ato+APQE+ ^PnPE^ (2) 

where to is the time of first maximum, Ato is the uncer- 
tainty in this time, Pq is the pulsation period at this time 
of first maximum and APq is the error in the observed 
period. Thus, any secular change in the period, dP/dt, 
will cause a parabolic curvature in an (O — C) diagram. 

We construct an initial (O — C) diagram using the pe- 
riods identified in Section IStI We then iteratively adjust 
to and Pq by the zeroth- and first-order terms from our 
best second-order fit until the adjustments are smaller 
than the uncertainty in these terms; these uncertainties 
result from the covariance matrix. Our recomputed, final 
(O — C) diagrams use these final ephemera and periods, 
which can be found for each mode in Table [H The fi- 
nal (O — C) diagrams are shown in Figure [2j Table [3] 
presents the times of maximum from each subset of our 
observations. 

3.3. Observed Rates of Period Change 

iMullallv etin (|20"ol ) published an (O - C) diagram 
of /i in WD 0111+0018, noting the large curvature but 
urging caution given the sparse coverage. Still, a rather 
high rate of change of period with time was observed: 
(3.87 + 0.43) X 10^12 s g-i 

We have added another 86 hr of observations, more 
than doubling the coverage, and confirm this large trend 
in the (O — C) diagram of /i. We have additionally been 
able to make (O — C) diagrams for three other periodic- 
ities. These diagrams can be found in Figure [5] 

Using Equation 2, a best-fit parabola to each (O — C) 
diagram yields a rate of change of period with time. We 
list these results in Table d Our results for /i a ctually 
fall within the errors of the lMullallv"eran (pOOl result, 
at (4.280+0.042) x IQ-^^ s s'^ (or 0.135+0.001 ms yr"!). 
We also strongly establish the need for a second-order fit 
rather than a first-order fit to the 18 data points: The 
best second-order fit for the {O — C) diagram of /i has 
— 145 (15 degrees of freedom, or d.o.f.), whereas the 
best linear fit has = 10292 (16 d.o.f.). 

Although the (O - C) diagrams for WD 0111+0018 
appear parabolic, as we would expect if the pulsations 
periods were changing secularly, we would also like to 
establish the timescale of variability if these observed 
trends are instead periodic in nature, since a portion of 
a sinusoid can mimic a parabola. For fi we find that the 
data can be well fit by a sinusoid with a 101 year period 
and a 3931 s (O - C) amplitude, with a = 142 (15 
d.o.f.). However, we will only continue our analysis in 
the framework of a secular period change. 

In the future, a third-order fit including the accelera- 
tion in the period change may be needed to best represent 
the data. This would enter the (O — C) equation as 

O - C = Ato + APaE+ ^PoPE^+ ^{PP^ + P^Po)E^ 

Currently, a third-order polynomial fit to the (O — C) 
diagram of /i has — 143 (14 d.o.f.), which reduces 
to a value still slightly higher than that for a second- 
order fit. The best third-order fit would correspond to 
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TABLE 2 

Mode periods and observed rates of period change in WD 0111+0018 



Mode 


to 


Po 


Amplitude 


dP/dt 




(BJDtdb) 


(s) 


(%) 


(10-12 s s-i) 


h 


2452666.562825(21) 


292.944305(9) 


2.573 


4.280(42) 


/2 


2452666.562314(31) 


255.663971(2) 


1.428 


0.340(58) 


2/l 


2452666.562701(42) 


146.472153(3) 


0.633 


2.075(44) 


h + h 


2452666.562624(29) 


136.518707(5) 


0.730 


1.109(30) 



600 F 



400 - 



200 - 



I/' -h 



TP~=\ 



u 

I 

o 



fj: 292.9 s 

dP/dt = 4.280(42)E-12 s s"' 
P/P = 2.169(21) Myr 



1 2fi: 146.5 s 

dP/dt = 2.075(44)E-12 s s" 
F/P = 2.237(46) Myr 



-I— I- 



I I I I I I I I I I I I I I I I I _ M 1 I I I I I I I I I I I I I I I I 



300 - 



^ 200 
o 

I 

o 

100 - 



l^: 255.7 s 

dP/dt = 0.340(5B)E-12 s s"' 
P/P = 23.8(3.5) Myr 



JL fi + f2: 136.5 s 

-- dP/dt = 1.109(30)E-12 s s"' 
-- P/P = 3.90(10) Myr 



-a ■ 

_i i_ 




_L 



_L 
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Fig. 2. — (O — C) diagrams for the four highest-amplitude periodicities present in WD 0111+0018. The /i and 2/i (O — C) diagrams 
are nearly identical, which strongly suggests that nonlinear combination frequencies in DAVs are not independent pulsation modes but are 
directly tied to their parent modes. A best-fit parabola yields a rate of change of period with time: See Table [2] for full solutions. 



with a second-order term 

12 „ o-l 



P= (-2.3±1.4)xl0-^iss- 

corresponding to dP/dt = (4.65 ± 0.23) x lO"^' s s 

We have also attempted a direct measurement of dP/ dt 
by computing a nonhnear least-squares fit to find the 
period for each subgroup. A best-fit line to these periods 
over time yields a rate of change of period with time. 
Unfortunately, there is simply not enough data for this 
direct measurement. The best result comes for /i, where 
we determine dP/dt = (5.5±6.8) x 10^^^ s s""'^, consistent 
with the (O — C) value but not yet significant. 

The other independent parent mode in 
WD 0111+0018, /2, shows a much slower period 
change than /i. Still, it is significantly faster than 
predicted by cooling theory, with a best-fit parabola 
that corresponds to (3.40 ± 0.58) x 10"^^ s s~^ (or 



0.011 ± 0.002 ms yr In this second-order fit 

is only marginally better than a linear fit: The best-fit 
parabola has = 93.8 (15 d.o.f.), whereas the best-fit 
straight line has = 128.7 (16 d.o.f.). 

Noise can cause scatter about a straight line in an {O — 
C) d iagram, especially over a relatively short baseline 
(e.g.. lKeDler et aLllTool . Since a bad period can cause a 
dominant linear term, we caution that the rate of change 
of the /2 periodicity may still be somewhat slower than 
the value cited in Table [5] if the period we use for our 
[O — C) diagrams is off by as little as 0.06 ms. Therefore, 
we recommend that our value be construed as an upper 
limit on the rate of change of period with time of /2. As 
with the other parent mode, we also investigated the best 
periodic rather than secular fit, and found a sinusoid with 
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TABLE 3 

Observed times of maximum in WD 0111- 



0018 



Time of Maximum 


Epoch 


(O-C) 


cr 


(BJDtdb) 




(s) 


(s) 


/l 2452668.431022 


551 


-0.11 


4.41 


2452968.919383 


89176 


5.29 


1.96 


2453353.385279 


202569 


25.11 


1.62 


2453709.441810 


307583 


56.14 


1.30 


2453725.336726 


312271 


53.99 


1.34 


2453999.274020 


393065 


94.00 


1.38 


2454064.667839 


412352 


103.15 


1.33 


2455067.170353 


708026 


305.96 


1.09 


2455089.880472 


714724 


319.29 


1.69 


2455152.022733 


733052 


327.43 


0.87 


2455209.164001 


749905 


342.55 


1.31 


2455452.386565 


821640 


412.39 


1.24 


2455477.944792 


829178 


429.02 


1.00 


2455511.887808 


839189 


440.13 


1.29 


2455565.211171 


854916 


443.63 


1.37 


2455837.226626 


935143 


536.16 


1.07 


2455857.749832 


941196 


549.35 


1.67 


2455959.724574 


971272 


574.11 


2.94 


/2 2452668.432414 


632 


-2.94 


6.75 


2452968.917516 


102179 


0.57 


3.07 


2453353.384106 


232107 


5.55 


2.53 


2453709.440530 


352434 


1.93 


2.15 


2453725.333732 


357805 


3.46 


2.00 


2453999.273126 


450381 


19.25 


2.58 


2454064.668563 


472481 


11.27 


2.06 


2455067.167424 


811269 


27.44 


1.75 


2455089.875468 


818943 


37.12 


2.11 


2455152.018898 


839944 


30.44 


1.34 


2455209.158663 


859254 


34.84 


2.40 


2455452.382727 


941450 


38.19 


1.96 


2455477.937402 


950086 


48.08 


1.66 


2455511.880763 


961557 


33.06 


1.86 


2455565.206336 


979578 


42.18 


2.23 


2455837.222237 


1071504 


49.81 


1.65 


2455857.743562 


1078439 


62.61 


2.60 


2455959.716124 


1112900 


55.88 


4.21 


2/i 2452668.432466 


1103 


-11.08 


10.53 


2452968.917468 


178351 


-3.06 


3.98 


2453353.385256 


405138 


33.60 


3.11 


2453709.441750 


615166 


61.32 


2.54 


2453725.336700 


624542 


62.18 


2.86 


2453999.273975 


786130 


100.45 


2.92 


2454064.669411 


824705 


102.78 


2.50 


2455067.170145 


1416052 


297.94 


2.62 


2455089.878699 


1429447 


322.50 


2.64 


2455152.022654 
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Fig. 3. — The amplitude evolution of WD 0111+0018 over nine 
years of observations. The dotted lines represent the weighted 
mean amplitude over the entire data set, shown in Table [J] While 
our data show that the amplitudes are not stable over the entire 
data set, they rule out a large-scale amplitude increase as the cause 
of the large trends in the (O — C) diagrams shown in Figure [2] 



a 9540 year period and 3051490 s (O— C) amphtude, with 
a = 93.8 (15 d.o.f.). (If such a trend were from an 
external companion, it would require a 10.8 Mq unseen 
companion at 401 AU.) 

A Monte Carlo simulation o f the u ncertainties using 
the software package Pcriod04 (jLenz fc Bregcr 2005) in- 
dicates that our formal least-squares uncertainties on the 
phases used to construct the (O — C) diagrams are un- 
derestimated by at most 10%. Thus, the scatter about 
our residuals (evident by the large Xred quoted above for 
the second-order fits, e.g. y^^^ = 9.7 for fx) may not be 
entirely the result of underestimated uncertainties. The 
deviation about the residuals is consistent with the ob- 
served scatter in other long-studied DAVs , such as Gil 7- 
B15A (see Figure 2 of lKepler et ani2005h . 

We have investigated the stability of the four highest- 
amplitude periodicities in WD 0111+0018 by plotting the 
amplitudes of each determined by a linear least-squares 
fit to each subgroup, using the periods in Table [21 This 
result is shown in Figure [31 It is meaningful since all 
observations were obtained and reduced in an identical 
manner. We have also marked the weighted mean am- 
plitude for each mode with dotted lines. 

There are considerable deviations about this mean 
compared to the formal uncertainties, such that only the 
periodicity at 2/i has a constant amplitude in a statis- 
tically meaningful way. While the amplitudes of the /i 
and /2 modes do not appear constant over our entire data 
set, the deviations from the mean do not appear coher- 
ent, and cannot explain the large-scale parabolic changes 
we observe in observed times-of-maxima. 

Finally, in line with the original intentions of this data 
set, we have investigated whether periodic changes in the 
times of maxima can explain the scatter observed about 
the residuals after subtracting the dominant parabolic 
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terms. This scatter might reveal the presence of an un- 
seen companion. However, a Fourier transform of all 
four [O — C) residuals, after subtracting the dominant 
parabolic terms identified in Figure [21 show no peaks 
above twice the mean noise level of 1.5 s. We can thus 
conservatively rule out the presence of any companion 
more massive than Jupiter at a current orbit between 
2 — 10 AU; the scatter is predominantly incoherent and 
cannot currently be attributed to an unseen companion. 

4. DISCUSSION 

4.1. Nonlinear Combination Frequencies 

Before focusing on the rapid rates of period change, we 
first draw attention to the observation that dP/ dt for the 
2/i and /i + /2 combination frequencies evolve exactly 
in lockstep with their parent modes. For example, we 
expect the 2/i combination frequency to have half the 
rate of change of period with time of its parent mode, 
since 2/i has half the period. We expect (2.140±0.042) x 
10-12 g g-i and in fact observe (2.075 ± 0.044) x lO^^^ g 
s-i, in excellent agreement. 

The /i + /2 combination frequency is a bit more com- 
plicated, given it is a sum of the parent frequencies and 
not an integer ratio. We expect its rate of change of 
period with time to be 

dt p2 ^^ + p2 

where Pi = fr\ P2 - f2-\ and P3 = (/i + f2)-^- We 
thus expect a rate of (1.027 ±0.099) x 10" s s^^ and in 
fact observe (1.109 ± 0.030) x lO'^^ s s'K As with 2/i, 
our observed rates match the predictions exactly within 
the errors. 

This lends convincing evidence that these nonlin- 
ear combination frequencies are not independent modes 
themselves, but rather artifacts of some nonlinear distor- 
tion in the outer layers of the star. They are thus directly 
tied to the parent modes, and evolve at exactly the same 
rate as these parent modes. This is an important result 
that convincingly demonstrates the nonlinear nature of 
these signals. 

In fact, an extensive theory explaining the creation of 
nonlinear combinat i on fre quencies in DAVs was intro- 
duced by iBrickhilll ()1992|) . He attributed the nonlin- 
earities to the changing thickness of the star's convec- 
tion zone, which acts to distort a linear input signal. 
Physically, local surface temperature variations lead to 
changes in the depth of the convection zone, which both 
absorb and release energy. Because this pro cess is so sen- 
sitive to temperature (cx T~^^ for a DAV; iMontgomervl 
1200^ . these variations in the depth of the convection 
zone lead to nonlinear effects in the observed light curve. 
An alternate explanation t hat may play a role fo r the 
hottest DAVs is given by iBrassard et al.l (|1995[ ) and 
[Fontaine fc BrassardI ()2008[ ). They invoke the nonlinear 
response of the flux observed in a passband to the tem- 
perature (i.e., the "T^" nonlinearity) . Regardless of the 
source of these nonlinearities, our result confirms that the 
phase evolution of these combination frequencies mirrors 
that of their parent modes. These combination frequen- 
cies are thus directly tied to their parent modes. 

The presence of nonlinear combination frequencies may 
also be useful in identifying the spherical degree (£) and 




I values m values 



Fig. 4. — A probability distribution ol£ and m values for the best 
7% of fits using the amplitude ratio of the nonlinear combination 
fr equencies to the par ent modes, a method similar to the one used 
in lYeates et al.l II2005I ). The amplitude ratios for the 292.9 s mode, 
/l, are best explained by an £ = 1, m = mode (top two panels). 
The amplitude ratios for the 255.7 s mode, /2, indicate this is an 
£ = 2 mode, although there is more ambiguity about the azimuthal 
order of this mode (bottom panels). 

azimuthal order (m) of the parent modes. In an inves- 
tigation of eight DAVs, Yeates et al. (2005;) found that 
the amplitudes of combination frequencies are too large 
to be created by the nonlinearity, and are probably 
due to the higher-order temperature dependence of the 
convection zone. They employed the formalism of iWul 
()2001| ) , who derived an analytical expression for the am- 
plitudes of the combination frequencies in the context 
of the Brickhill framework of a depth- varying convection 
zone. Using the ratio of the amplitudes of the parent and 
combination frequencies, Yeates et al. derived constraints 
on the £ and m value of the parent mode. 

With the same goal in mind, we have used our own 
code to automatically search for combinations of £ and m 
that best fit the observed amplitudes in WD 01114-0018. 
This code uses a genetic alg orithm to efficien t ly sam ple 
parameter space; see §4.1 of [Provencal et all ()2012D for 
more details. We run this algorithm 1000 times so that 
we obtain 1000 e stimates o f the best-fit parameters. 

The method of lWul (|2001l ) is only applicable to second- 
order combinations, i.e., 2/i or /1 + /2 but not 3/i. Even 
with just three independent combination amplitudes we 
are able to find constraints. Consistently, we best repro- 
duce the observed amplitude ratios (we have used the am- 
plitudes of the 2009-2012 group, provided in Section |3T]) 
when we allow the /i parent mode at 292.9 s to be an 
^ = 1, m = mode. We also find the /2 parent mode at 
255.7 s to be an £ = 2, m = 0, ±1 mode. This is shown in 
Figure m where we plot the distribution of the top 7% of 
fits with the lowest x^j interpreting this as a probability 
distribution for the identifications. 

In addition, we have implemented the convcctivc non- 
linear light curve fi t ting c ode of , Mont gome ry (2005) and 
iMontgomerv et all (|2010( ) to generate a synthetic light 
curve based on the parent modes in WD 0111-1-0018. 
This method investigates the time-averaged thermal re- 
sponse time of the convection zone, Tc, which relates to 
the mass and depth of the convection zone. In princi- 
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pie, it sho uld be superior to the perturbative method of 
iWul (|2001f ) in that it makes use of all higher-order am- 
plitudes and their phases. These nonlinear fits produce 
estimates for the t and m values of the parent modes, 
the time-averaged convective timescale tq, the tempera- 
ture sensitivity exponent N , and the inclination of the 
pulsation axis to our line of sight 6. 

Fitting our light curves from the 2011b subgroup, we 
find an identical match to our previous method. The 
292.9 s mode is best represented by an ^ = 1, m = 
mode and the /2 parent mode at 255.7 s best fits as 
an ^ = 2, m = ±1 mode. Additionally, these fits find 
To = 137.8±11.8s, N = -82.0±3.6,and6' = 34.6±2.5°. 

We thus adopt an interpretation that /i is an m = 
dipole {£ = 1) mode and /2 is a quadropole {£ = 2) mode. 
The constraints on the azimuthal order of /2 are less 
conclusive. We will use this information to help constrain 
our asteroseismic models, which will in turn inform our 
attempts at understanding the observed dP/dt. 

4.2. Expected Rates of Period Change 

We can immediately rule out the possibility that an un- 
seen external body is vigorously tugging WD 0111+0018 
away from our line of sight, causing a light-travel time 
effect that is manifest as the increasing delay in pulse ar- 
rival times: Rapid space motion would affect all modes 
identically, and the periodic modulation in the mode /2 
is more than an order of magnitude slower than that of 
/i, as we found in Section [5751 We must be observing a 
phenomenon internal to the WD. 

We expect cooling to be the dominant evolutionary ef- 
fect of an isolated DAV, and this cooling timescale should 
have a direct effect on the period evolution of modes in 
a pulsating WD. To first order, as a WD cools it be- 
comes more degenerate. This causes a decrease in the 
square of the Brunt- Vaisala frequency, N"^. The Brunt- 
Vaisala frequency is essentially the oscillatory frequency 
of an adiabatic fluid displacement in a convectively stable 
medium. In the limit of low temperature and complete 
degeneracy, where pressure does not depend on tempera- 
ture, a displaced fluid element would remain in pressure 
and density equilibrium with its surroundings. It would 
stay put and thus have zero oscillatory motion — and 
thus a zero 7V^. Thus, decreases with increasin g de- 
gener acy, and therefore decreasing temperature fsee lKimI 
12001 . 

The period of a ^-mode pulsation i s inversely propor - 
tional to the Brunt- Vaisala frequency (jUnno et al.lll989l) . 
Therefore, cooling will cause a secular increase in a pul- 
sation period with time, which can be estimated from the 
cooling timescale predicted from Mestel theory, as well 
as full evolutionary models. Predicted rates of period 
change for DAVs undergoing simple cooling are of order 
10^^^ s s^^; we flnd observational evidence to support 
this timescale from observations of G117-B 15A and R548 
(jKepler et al.|[20Tlt [M^dam et al.l [20121 ). The 215.2 s 
mode in G117-B15A and the 213.1 s mode in R548 are 
most likely either an ^ = 1, fc = 1 mode or a i = 1, k = 2 
mode, depending on the adopted hydrogen l ayer mass 
(|Bischoff-Kim et al.l[200l [Romero et all [2012'). 

We note that the work of monitoring long-term phase 
evolution of DAVs requires a mode coherent in ampli- 
tude and phase. This is usually not the case for the 
cooler DAVs, which evidence strong amplitude modula- 



tion and a corresponding lack of stabili ty in phase (e.g. 
[Kleinman et"all [19981: [Polez et al.|[2006l) . One cool DAV 
has, for at least a brief time, shown enough phase co- 
herence to create a mea ningful (O — C) diagram: The 
615.2 s mode in G29-38 (Kleinman 1995) , which showed 
a coherent but extremely rapid phase change, more than 
200 s over less than three months of monitoring in 1988 
(jWinget et al.l [1990fl . There is also evidence that the 
274.8 s mode in R54 8 has a dP/dt > IQ-^^ s s'^ 
(|Mukadam et al.|[2003[ ). and that the 270.5 s and 304.1 
modes in G117-B1 5A have dP/dt > 3 x 10"" s s^^ 
(jKepler et al.|[2005l ). Still, most modes in DAVs observed 
to date simply do not show enough coherence to create 
stable (O — C) diagrams. 

WD stars with heavier cores should cool faster than 
those that are lighter. Indeed, spectral flts show that 
WD 0111+0018, at ~ 0.71 Mq, is shghtly more massive 
than G117-B15A and R548, which are between 0.59- 
0.60 Mq. StiU, even a 1.4 Mq Fe-core WD would have 
a dP/dt from cooling of less tha n 9 x 10~^^ s s~^ if the 
WD has Toff' = 11800 K (Bradlc fiFaI|[T99^ . 

Our observed rate of change of period with time for 
/i in WD 0111+0018 exceeds lO'^^ g g-i^ j^ore than 
two orders of magnitude faster than predicted from 
cooling alone. It is also more than an order of mag- 
nitude faster than expec ted from an avoided crossing 
(jBradlev fc Winge3[1991|) or by invoking excess cooling 
from exotic particles such as axions (Bisclioff-Ki m et al.l 
[20081: [Corsico et all [20121 ). We thus interpret this dis- 
covery as evidence of an additional physical timescale 
acting on the evolution of the pulsation periods in 
WD 0111+0018. 

4.3. Asteroseismology of WD 0111+0018 

Both the 292.9 s and 255.7 s modes of WD 0111+0018 
likely probe different regions of the star, so mode identi- 
flcation (specifically the radial order, k) would yield fur- 
ther insight into understanding these anomalous rates of 
period change. Unfortunately, asteroseismology is made 
more difficult by the presence of just two independent 
modes; there are far more free parameters than observa- 
tions. 

Still, we have attempted to match the observed peri- 
ods to adiabatic pulsation models with some assumptions 
and the constraints provided by our spectroscopic mass 
and temperature determinations. Our two best matches 
to the observed periods have the following properties: 

Solution 1: /i is an ^ = 1, fc = 5 mode (m — occurs 
at 293.4 s), while /2 is an £ = 2, fc = 8 mode (m = 
occurs at 255.0 s). This WD model has a stellar mass of 
0.710 Mq, an effective temperature Tos = 11630 K, a He 
envelope mass of A/ho/^^* = 10~^'^, and a H envelope 
mass of Mh/A/, = 10"*-^ 

Solution 2: /i is an £ = 2, fc 10 mode (m = occurs 
at 293.7 s), while /2 is an £ = 2, fc = 8 mode (m = 
occurs at 255.0 s). This WD model has a stellar mass of 
0.720 Mq, an effective temperature Tcff = 11610 K, a He 
envelope mass of Mho/M* = 10"2 °, and a H envelope 
mass of A/h/M* = lO'^'-^ 

Additionally, we can u se our spectroscopic constraints 
to explore the models of [Romero et al.[ (j2012D of a 0.705 
Mq, 11810 K WD with a thick (10-4-445 m^) hydrogen 
layer mass. Similarly, we find that the 292.9 s mode is 



8 



Hermes et al. 









p-modes 






- 
- 




2 


- 








- 






- 










- 







- 










- 
































- 


















o 


-2 


- 










- 








/ 250 sec 












-4 




300 sec 
















g-modes 












-6 




, , , 






1 , 





-5 -10 -15 

log(l-MyMj 



Fig. 5. — A propagation diagram for a representative asteroseis- 
mic fit, a WD model with stellar mass 0.710 Mq, an effective tem- 
perature Toff = 11630 K, a He envelope mass of Mue/A^* = lO"'^ '^, 
and a H envelope mass of Mjj/M* = 10~*'^. The run of the Brunt- 
Vaisala frequency is shown in blue. The run of the Lamb (acoustic) 
frequency for an £ = 1 mode is shown as a dashed magenta line, 
while the Lamb frequency for an £ = 2 mode is shown as a solid 
black line. The horizontal axis is in fractional mass units; the sur- 
face is to the far right, while the center is to the far left. The base 
of the convection zone can be seen as the sharp drop in A'^^ around 
(1 - A/r/A/,) ~ 10-1''. 

either an ^ = 2, fc = 10 or £ = 1, fc = 5 mode. The 255.7 
s mode is either an ^ = 2, fc = 8 or ^ = 1, fc = 4 mode. 
We caution that because these solutions are consistent 
does not imply they are correct — our assumptions may 
not indeed be valid. But these solutions will help inform 
our analysis of what could be causing the rapid period 
change in WD Olll-f-OOlS. 

Using our spherical degree {€) identifications described 
in Section 14.11 we can break the degeneracy in /i and 
suggest that the 292.9 s periodicity represents an £ = 
l,fc = 5 mode, while the the 255.7 s periodicity is an 
£ = 2, fc = 8 mode. In order to visualize how each such 
mode samples the star, we have included a propagation 
diagram in Figure [5] We have marked where a ~250 s 
and a ^300 s ij-mode for different t values would probe 
the star (where cr^ < A^^, L\). 

However, an even more informative display of how 
each mode could probe the star differently can be 
found in Figure [6l Here we plot, a s a function of the 
"normalized buoyancy radius" (see 'Montg omerv et all 
20031). the weight fu nctions for A^^ (see Equation 8c of 
Kawaler et all [19851) for our best identifications of the 
modes in WD 0111-1-0018. The weight function essen- 
tially shows where in the star an excited pulsation is most 
resonant; it shows the sensitivity of a mode's period to 
changes in the equilibrium model. 

Immediately evident from this exercise is that the ma- 
jority of the energy of the 292.9 s mode, /i, is located in 
the He layer or the region of the He-H transition. Con- 
versely, nearly all the energy of the 255.7 s mode, /2, 
is resonant much deeper in the star, inside the C/O-He 
transition. This suggests that the 255.7 s mode energy is 
predominantly resonant inside the core, while the 292.9 
s mode samples significantly more of the surface layers 
of the star. The 255.7 s mode would therefore be much 
less sensitive to changes in the outer regions of the WD. 
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Fig. 6. — We can investigate how each pulsation samples the 
star by plotting a weight function (WF$) for each mode. The top 
panel shows the chemical transition zones as a function of depth in 
the star. The bottom two panels show the weight functions for our 
best models (which are only marginally constrained) to represent /i 
and /2 , respectively. All plots are shown as a function of th e "nor- 
malized buoyancy radius," "l?(r), described in ^Montgomer y et al.l 
(2003), where is the center and 1 is the surface of the WD. The 
corresponding mass fraction is shown at the top of this figure. 

4.4. Additional Mechanisms for Rapid Period Change 

Given the inscnsitivity of the 255.7 s mode to changes 
in the outer layer of the star, one possibility for the new 
timescale represented by the fast rate of change of pe- 
riod with time in WD 0111-1-0018 involves a change in 
the surface layers of the star at a rate much faster than 
from cooling alone but much slower than the timescales 
associated with the pulsations themselves. For instance, 
the depth of the convection zone of the WD could change 
gradually, perhaps as a result of a change in the global 
magnetic field of the WD. We have seen evidence for 
rapid changes in the pulsation profile of WDs in the 
past: especially notable is the so-called "sforzando" event 
seen in GD 358 in 1996, in which this DBV changed 
from a highly multi-periodic pulsation spectrum to a 
higher-amplitude, monoperiodic pulse shape, indicating 
a rapid thinning of the convectio n zone during t his event 
(Provencal et al. 2009; Montgom erv et al.l[2010[ ). 

We have investigated the plausibility of this scenario by 
comparing two slightly different theoretical models. Both 
model WDs have identical masses (0.71 M©), effective 
temperatures (11810 K), and chemical profiles. However, 
we have changed the efficiency of convection by varying 
the mixing length parameter: in our initial model we use 
ML2/a=1.23 and in the other we use ML2/a=1.12. This 
effectively, albeit superficially, causes the second model 
to have a slightly thinner convection zone. 

We take the difference of the periods computed by the 
first model with the periods from the second to estimate 
the effect these changes would have on different modes. 
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We find that even the shortest period I = 1, A: = 1 pul- 
sations would change in mode period by more than 1 ms 
from one model to the next. This result can reproduce 
the observed rates of period change in WD 0111+0018 if 
the timescale of this change in the depth of the convec- 
tion zone acts on the order a few decades. 

One way to perturb the overall depth of the convec- 
tion zone is to introduce a gradual change in the global 
magnetic field strength. In fact, increasing the global 
magnetic field from G to 1 kG would effectively change 
the depth of the convection zone by the same amount 
as changing the convective efficiency from ML2/a=1.23 
to ML2/a=1.12 (Montgomery & Vishniac, in prepara- 
tion). It is possible then that the rapid period change 
we observe in WD 0111+0018 could be related to a secu- 
lar (or long-term periodic) change in the strength of the 
magnetic field in the star. If this is indeed a byproduct 
of a stellar cycle, we might expect our heretofore secular 
increase to at some point turn over or in fact become dis- 
continuous when the magnetic cycle changes course. Fu- 
ture observations of WD 0111+0018 will help constrain 
this hypothesis. 

Another mechanism that could affect one mode differ- 
ently than another, in the way we see in WD 0111+0018, 
is a changing rotation rate. Rotation will act to break 
spherical symmetry, splitting non-radial pulsations into 
{21 + 1) azimuthal components that would norm ally be 
degen erate with the m = component (Hanse n et al.l 
119771) . However, if the rotation rate were somehow to 
change (a spin- up or spin-down) , only the m 7^ compo- 
nents would be affected. If /i is an m = ±1 mode and ji 
is an m = mode, a changing rotation rate could explain 
the large difference in rates of period change between the 
two. 

Our observed rate of change of period with time for 
/i would require the rotation of the WD to change at 
a rate |Prot| > 0.3 s yr~^. Such a rapid change would 
require a tremendous reservoir of angular momentum to 
redistribute, and some place to deposit this energy. This 
appears especially unlikely given the observational result 
that the hot pulsating white dwarf PG1159-035 rotates 
as a solid body through more than 97.5% of its mass, 
indicating that WDs hav e very little angular mo mentum, 
even at T^s > 80,000 K (jCharpinet et al.|[2009l ). 

Unfortunately, we do not observe any other multiplet 
components for fi or /2 , making it difficult to definitively 
rule out this rotational effec10. But our results from the 
convective light curve fitting discussed in Section |4?T] did 
yield predictions for the m value of each mode: m — 
for /i and m = 0, ±1 for /2. This contradicts the values 
needed to invoke a changing rotation rate to explain the 
different rates of period change between these modes. 

Some additional effects intrinsic to the star that could 
cause phase variations include linear interference (such 
as beating of closely spaced modes), nonlinear interac- 
tions with other modes (such as mode coupling) and 
a precession of the rotation axis with our line of sight 
(jWinget et al.lll99"ot ). Each of these cases would simulta- 

^ Due to inclination effects, it is not uncommon to see an m = ±1 
component of a rotationally split multiplet and not the m = 
central compo nent. This is well illus trated for different £,m values 
in Figure 1 of IBrassard et al] 1199511 . It is also possible that only 
one component of a multiplet is excited. 



neously affect the amplitude of a given mode. Since we 
do not see evidence of the amplitudes evolving identically 
with the phase (comparing Figure [2] with Figure [S]), none 
of these are likely responsible for the inferred high rates 
of period change in WD 0111+0018. 

5. CONCLUSIONS 

As part of a search for substellar companions to iso- 
lated WDs, we have discovered a DAV, WD 0111+0018, 
with an unexpectedly fast rate of change of period with 
time. We see this anomalously high rate reflected in not 
only the highest- amplitude mode in the star at 292.9 s 
but also the two nonlinear combination frequencies that 
sample this parent mode. All three periodicities are 
changing at a rate faster than 10~^^ s s~^, more than two 
orders of magnitude more quickly than expected from 
cooling alone in a DAV. 

The pulsation period of the other parent mode at 255.7 
s is increasing at a rate at least an order of magnitude 
more slowly, which allows us to rule out an extrinsic 
source of this rapid period change. Instead we are likely 
observing a physical phenomena intrinsic to the star act- 
ing on a timescale we have not yet observed in a DAV. 

This discovery indicates that our understanding of the 
expected rates of period change for the slow and simply 
evolving DAVs is far from complete. It also complicates, 
in a macroscopic sense, our ability to infer the contribu- 
tions of exotic weakly interacting particles such as axions 
to anomalies in the rate of change of period with time of 
DAVs. 

We have attempted to use the nonlinear combination 
frequencies present in this star as probes of the spherical 
degree of the two parent modes . Using convective light 
curve fitting (|Montgomervll2005D . we find that the 292.9 
s mode is best represented as an £ = 1, m = mode and 
that the /2 parent mode at 255.7 s best fits as an € = 2, 
m = ±1 mode. This identification, if confirmed, would 
rule out a changing rotation rate as the explanation of 
the rapid period change in the 292.9 s mode. 

Guided by the spherical degrees we identify by this con- 
vective light curve fitting we have found a best-fit aster- 
oseismic model, which is admittedly poorly constrained. 
Still, this exercise suggests that fi is an £ — l,k — 5 
mode and /2 is an £ = 2, fc = 8 mode. We construct a 
weight function, and show that this solution for /2 is al- 
most entirely resonant inside the core, and is thus much 
less sensitive than /i to changes in the outer regions of 
the WD. 

We discuss a scenario where a gradual thinning of the 
convection zone could be responsible for the rapid rates 
of period change observed in this DAV. This thinning 
could be caused by increasing the global magnetic field 
by less than 1 kG over a decade. Such a hypothesis is 
hard to test empirically, although if this magnetic field 
change is related to a stellar magnetic activity cycle, we 
would expect an eventual turnover or discontinuity in the 
observed (O — C) diagram. 

Independent of the cause of this unexpected period 
change, we have for the first time shown that the evo- 
lution of nonlinear combination frequencies in a DAV 
matches the evolution of their parent modes. This behav- 
ior confirms that these combination frequencies are not 
independent pulsation modes but rather artifacts from 
some nonlinear distortion occurring in the star. 
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